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Abstract 



This research is concerned with the stability and time-dependent 
motion of a spinning satellite, simulated by a rigid body with a 
cavity partially filled with liquid. The work includes the problem 
formulation, consisting of the boundary-value problem for the liquid 
and moment equations for the entire system. Because of large Reynold's 
numbers involved, viscosity effects are negligible everywhere except 
for a thin boundary layer near the wetted surface. Using a boundary- 
layer analysis, the effect of the boundary layer is replaced by modified 
boundary conditions for the liquid. The solution of the differential 
equations for the inviscid problem has been solved in closed form. A 
semi -analytical numerical solution of the inviscid equations subject to 
the viscous boundary condition has proved unsuccessful. 
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1 . Kineiiiatica l Relations 

Let XYZ be an Inertial coordinate system with the origin at the 
mass center c of the satellite (see Fig. 1), and let the satellite 
be rotating with an angular velocity u with respect to axes XYZ. 
Oenoti.ig by xyz a set of body axes, and by T, J, Ic the unit vectors 
along these axes, the angular velocity of the satellite can be written 
in the form 


u) ® 0) T + 10 T + (0 7 ( 1 ) 

A y z 

where (o , lo , (o are the corresponding angular velocity components. 

A y z 

It will prove convenient to introduce a set of axes Cnc parallel 
to axes xyz but with the origin at the center 0 of the spherical tank 
instead of the center c of the satellite, as shown in Fig. 1. Denoting 
by Rq the vector from c to 0 and by r the vector from 0 to any position 
in the li'iuid, the absolute position of any point in the fluid can be 
written as 


R = Rq + r = xi + yj + zk 

where 


( 2 ) 


«0 ' * ''o/J * »> 

and 

^ + nj + Ck (4) 

Differentiating Eq. (2) with respect to time, and recognizing 
the fact that the unit vectors T, J, k rotate with angular velocity w. 


1 


1 


we obtain the absolute velocity at <• point In the liquid In the form 
where 7 q Is the velocity of the origin 0 and 

( 6 ) 

is the velocity of the liquid relative to the tank. The absolute 
acceleration of the liquid can be obtained by differentiating Eq. (5) 
with respect to time. The result Is 


a 


dv 

at 




0) X q + 


clw 

at 


X r + 0 ) X 


dr 

at 


(7) 


where 


Dt 



V 


( 8 ) 


Is known as the substantial derivative operator. Denoting dw/dt by 
w, and inserting Eq. (8) Into Eq. (7), we can express the absolute 
acceleration at a point In the liquid as 


a = ao + |^+q'Vq + 2u)xq + a)Xr + a)x(u)Xr) (9) 

in which ^ Is known as the convective acceleration. The absolute 
acceleration of a point in the rigid part of the satellite can be obtained 
by letting ^ = 0 in Eq. (9). 
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2. The Mathematical Formuidtion for the Motion of the Sy»tem 

The problem formulation consists of the boundary-value problem, 
namely, the differential equation and boundary condition for the liquid, 
and the moment equation governing the rotation of the satellite. 

2.1 The equations governing the liquid motion 

For an Incompressible fluid, the velocity vector q must satisfy 
the continuity equation 


V • q ■ 0 (10) 

In this case, the Navler-Stokes equations become 

a‘«?Q + |^+q • Vq + 2wxq + o)Xr + a)’x (wxf) 

- - 1 yp + 7 + V V^q (11) 

where 7 Is the body force per unit mass, and p, p, and v are th- 
density, pressure, and kinematic viscosity, respectively. 

From vector algebra, we have 

q • Vq = ^ V(q • ^ - q X (7 X q) 
u) X (to X r) « - ^ 7[(w X r) • (to x r)] (12) 

Iq = ^(^0 • r) 
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where V x q is known as the vorticity of the liquid relative to the 
tank. Assuming that the body force is derivable from a potential 
function, we can write 

7 ■ - VB .13) 

where B is a scalar function. It is convenient to introduce a pressure 
function p* defined by 

p* . i p 4 B + a^ r 4 1 q • q - ^ (w X r ) . (w A r ) (14) 

Then, Eq. (11) can be rewritten as 

l^-^x (Vxq) + 2wxq + ajxF» - Vp* + vV^q (15) 

2.2 Boundary conditions 

The boundary conditions on the liquid are of two types, namely, 
geometric and natural. The geometric boundary conditions are the vesult 
of kinematic considerations, whereas the natural boundary conditions are 
the result of dynamic considerations. Moreover, we must distinguish 
between the boundary conditions on the wetted surface and the free 
surface S^. 

At the wetted surface the velocity of the liquid relative to the 
satellite is zero. If W, t^ and I 2 denote unit vectors normal and 
tangential to the wetted surface, then we have 
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q ■ Q • n ■ 0 on S. (16a) 

n w 

‘^tl " ^ ‘’tZ • q • ^2 “ ° 

The boundary condition (16a) Is a result of the Impermeability of the 
rigid part of the satellite and must be satisfied whether the fluid Is 
viscous or Inviscld. On the other hand, the boundary conditions (16b) 
are referred to as the no-sllp conditions because. Irrespective of how 
small the viscosity of the fluid Is, there can be no relative tangential 
velocity between the fluid and the wetted surface. If the fluid Is 
Inviscld, the boundary conditions (16b) may not be satisfied. 

At the free surface there are three boundary conditions. The 
first Is kinematical In nature, expressing the fact that a particle on 
the free surface must always remain on the free surface that Is, If 
F(x,y,z,t) « 0 Is the equation of the free surface, the we must have 

^ [F(x,y,z,t)] » 0 on F(x,y,z,t) ■ 0 (17) 

Recalling definition (8), Eq. (17) becomes 

It * "x I; * V ly * "z ‘ 

where q , q , q are the components of ^ along the xyz axes. The other 
X y z 

boundary conoltiuns are dynamic, representing the continuity of normal 
and tangential stresses across the free Interface. Taking the surface 
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tension into account, we can write the continuity of the normal 
stresses 


T(K^ + K 2 ) + P ■ o on F(x,y,z,t) ■ 0 (19) 

where T is the surface tension, and K 2 are the curvatures of the 
free surface along two orthogonal directions, P represents the pressure 
in the vacant part of the tank, and o is the normal stress. The 
continuity of the tangential stresses is expressed by 

* 0 , T 2 ■ 0 on F(x,y,z,t) (20) 

2.3 The moment equation for the satellite 

We assume that the mass center c of the satellite is fixed in 
an inertial space, and, moreover, that it does not shift with respect 
to the rigid part of the satellite, even though the liquid may be 
moving. Recognizing that the acceleration of any particle in the 
rigid ,iart of the satellite can be obtained by letting ^ ■ 0 in Eq. 

(9), we can write the moment equation about c in the form 


R’x[ao + |^+q -Vq + Zwxq + uxr + wx (wx r)]dm^ + 


m. 


m.. 


R X [aQ + o)xr + wx (tox r)]dm^ 


) 


I R X [7(1 p B) + vV^^dm^ + 


m 


m 


R X T dm 
r r 


( 21 ) 
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where the subscripts t and r d?s1gnate the liquid and rigid parts of 
the satellite, respectively. Because the acceleration of the point 
0 Is 

a Q ■ w X + w X (w X ^q) (22) 
and recalling that + 7 ■ IT, we can write* 


Rx[u)xl^+ujx(a)x IT)Jdm- + R^x[u)xl? + u)x(u)x I^)Jdm 

t 

■T*u) + u)xT"u (23) 

where T Is the Inertia dyadic of the satellite as If It were entirely 
rigid. Using Eqs. (12), (14), and (23), we can rewrite Eq. (21) as 


T*u) + u)xT'U)»-| Rx {v[p* - a^ • r + ^ (u) X r) • (u) X r)] - 

qx (7xq)+|^+2o)xq + v7^q)dm^ + j x dm^ (24) 

•"r 

where the first Integral on the right side of Eq. (24) Is the disturbing 
torque due to the liquid motion and the second Integral Is due to body 
forces acting on the rigid part. 

3 Steady-State Equilibrium Configuratio n 

First we seek an equilibrium configuration In which the satellite 
rotates with a constant aniular velocity n about the Inertial axis Z 


* See Ref. 1 , Sec. 12.8 
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and In which the liquid Is at rest with respect to the satellite. More- 
over, the body force Is assumed to be zero. To describe this configu- 
ration, we Introduce a coordinate system x'y'z* rotating with the 
angular velocity n about the axis Z, and denote by T' , J', t* the unit 
vectors along these axes. Hence, the equilibrium configuration Is 
given by 

q-0, w-nF', U)»0, FQ(x,y,z) ■ 0 (25) 

where FQ(x,y,z) ■ 0 Is tht equat'oi of the free sur "e. 

Substituting Eq. (25) Into Eqs. (14) and (15), and recalling 
Eq. (22), we obtain 

P* “ p P' ~ \ y'^) “ const (26) 

The boundary conditions (16), (18), and (20) are automatically satisfied, 
whereas, boundary condition (19), together with Eq. (26), yields 

T(K^ + K^) + ^ P n^(x'^ + y'^) - const (27) 

which describes the free surface Fq(x', y', z') ■ 0. Hence, the free 
surface Is the Intersection of the curve given by Eq. (27) and the 
sphere 

(x' - t (y' - + ( 2 ‘ - • Rj (28) 

where Is the radius of the spherical tank. In the case In which the 
surface tension T Is negligible, the free surface Is given by the Inter- 
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section of the cylinder 

x'2 + y'2 - rJ (28) 

and the spherical tank given by Eq. (28), where the radius R of the 
cylinder Is a function of the volume of liquid In the tank. 

4. Linear Perturbation Problem 

Next, let us consider the case of a small perturbation from the 
steady-state solution discussed In Section 3. To this end, we write 


U) ■ 

nk' + cu)i + ... , u) “ + ... 

(30) 

- 

R' + e 

+ ... 

(31) 

q ■ 

eq^ + 

• • • 

(32) 

p " 

p' + e 

+ ... 

(33) 

F(x 

.y.z,t) 

■ R - R^ + e F-|((j),z,t) + ...» 0 

(34) 


where c Is a small parameter and R, ()), z a cylindrical coordinate system 
with the origin at c. Recognizing that the perturbed body Is Inclined 
slightly with respect to axes x'y'z', we let the body axes triad xyz be 
obtained from x'y'z' hy means of the small rotations c0^ about x', c 02 
about y', and £0^ about z', which enables us to write 

X * x' + e (o^y* - 02Z') + ... 

y ■ y' + e (0^z' - 02 x') + ... (35) 

z ■ z' + e ( 02 ^' - 02 y') + ... 
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so that, from Eqs. (31) and (35), It follows that 

^ ^ ^ 

Substituting Eqs. (30) - (34) into Eqs. (14), (15), (19) - (20), 
and (24), and equating coefficients of equal powers of c, we obtain 
the steady-state solution of Section 3, corresponding to the 0(e°) 
problem, and the following equations for the 0(e) problem 

V • q ^ - 0 (37) 

3qi _ ^ 

^ + 2J2l( X q.j + X R"' » 

-V[^ P] ■ ^ ^ X R^^)J + V V^q^ (38) 


T • u)^ + ftlc' X T • u).j + 0)^ X T • nF' » 


[R- X (1 Vp^ + ^ + 2fik X q^ - V V^q^) + R] x i Vp’]dm^ 


m, 


where q^ Is subject to 

the following boundary conditions 

(39) 

q^ « 0 

on 

r-«s 

(40a) 

% " ■'b “ “ 

9F, 

'’r * ° 

on 

r ■ Rj 

(40b) 

'on 

R - R^ 

(41) 
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"Rz 

^R4> 


Pi ^ ^ ° 


'’“‘sT * 3T> ■ “ 


‘>- I ’•'r. 


k <R^> * B sT^ 


0 


on R ■ R 

c 

on R ■ R 

c 

on R ■ R„ 
c 


(42a) 

(42b) 

(42c) 


where r, a, B is a spherical coordinate system centered at 0, whereas 
R, 4>, z is a cylindrical coordinate system centered at c. Note that 
the boundary condition (41) was obtained by inserting £q. (34) into the 
kinematic condition (18), performing a so-called "t'^ansfer of boundary 
conditions", and keeping the 0(e) (rms*. 


5. Method of Solutio n 

5.1 The equationi of motion in dimensionless form 

To solve Fqs. (37) - (39) subject to the boundary conditions 
(40) - (42), we first consider the fluid problem. To this end, we 
introduce the dimensionless quantities 


X *= x'/2R^ , y = y'/2Rj 


t — t n “ q / u)|^ 

A A 

fl ' 


q^ = q^ U e 


A A 

xt 




* See Ref. 2, Sec. 2.1.3 


z - z'/2Rj 


(43) 
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“ 0^ 1 + 02 7 + 63 k * Aojw 0 e 


" f -At 
N 


■ i !!- e»‘ 

1 1 

where U ■ 2R^u)j^ is a fictitious velocity and * Jl(C-A)/A ■»$ the 
nutational frequency of the rigid satellite, in which A is the principal 
moment of inertia about x or y and C the principal momeiit of inertia 

A 

about z. Moreover, A is the dimensioiilcss eigenvalue sought. Sub- 
stituting Eqs. (43) into Eqs. (37) and (38), we obtain 


V • q, = 0 


(44) 


AA A A AAA ^ 

Aq^ + 2fik' X q^ + A^^0 x R' = 


A A A A A 


- v[p^ - (riF X R')-(A0 X R' + file' X R^)] + ^V^q^ 

e 

where, in terms of cartesian components. 


(45) 


5 =T^+j^+l^^ 

8x 9y 3z 


represents a dimensionless nabla, and R^ = 4R^ is Reynold's 


number. 


In a similar fashion, the moment equation, Eq. (39), becomes 


A A A A 


A A A ^ A A ^ 

A T • 0 + ACl(lc ' xT-0 + 0xT-lc') = 


f ^ A A A A A A 1 A OA A A A ^ 

1 [R* X (Vp^ + Aq^ + 2nF / ^ 

Jm^ e 


(46) 
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— — 2 

where I = I/4R^ni Is the dimensionless Inertia dyadic of the complete 
sybtem as if the system were entirely rigid, in which m is the total 
mass (rigid plus liquid) of the system. Moreover, iri^ = m^/m is the 
dimensionless mass of the liquid. For simplicity, T is evaluated by 
regarding the mass of the liquid as frozen in the equilibrium con- 
figuration. 

Typical values of R^, and v show that R^ = O(IO^) or larger*. 
Hence, one can neglect the last term in Eq. (45) and obtain 

/S /S /\ /V AO/S a a a a a a a a a ^ 

Aq^ + ZnV X q^ + X 0 x R' = -V[p^ - (file x R ' ) • (X0 x + flic' x R^)] 

(47) 

However, Eq. (47) is of first order in the spatial variables rather than 
second order. Hence, the general solutions of Eqs. (44) and (47) cannot 
be expected to satisfy all the boundary conditions (40) - (42). Because 
letting -> « is equivalent to assuming an inviscid fluid, the no-slip 
boundary conditions (40b) and the shear boundary conditions (42b) and 
(42c) cannot be satisfied in general. 

To obtain a solution valid everywhere, we must supplement the 
solutions of Eqs. (44) and (47) by two boundary layers, one near the 
free surface and the other near the wetted surface, by using the method 
of composite expansions.** 

The result of the boundary layer analysis will be the modification 
of the boundary conditions for the inviscid problem by accounting for 


* See Ref. 3 

** See Ref. 2, Sec. 4.2 
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the liquid viscosity. Hence, the problem will reduce to the solution 
of an inviscid problem, with the effect of the viscosity reflected 
in the modified boundary conditions. 

The nature of the geometry of the liquid boundaries demands the 
use of a cylindrical coordinate system with the origin at c to de- 
scribe the boundary layer near the free surface and a spherical coordi- 
nate system with the origin at 0 to describe the boundary layer near 
the wetted surface. Note that either coordinate system can be used 
away from the boundariei . As a result, we shall seek the solution in 
two parts, one valid everywhere except near the free surface and the 
other valid everywhere except near the wetted surface. 

5.2 The boundary layer next to the free surface 

We assume the solution next to the frea su> face to have the 

form 

q^ = qj^e^ + q^e^ + q^e^ = [u^(R,(j>,z) + 6u^(o,()),z)]e,^ + 

[v^(R,(j),z) + v^(3,(j>,z)]e^ + [w.(R,4.,z) + Wy(p,<(),z)]e^ (48) 

Pi = Pi(R,<|).z) + Py(p,4>.z) (49) 

A A A 

where p = (R - R_)/6, in which 6 = l/»^ is proportional to the boundary 
layer thickness. Note that the subscripts i and v refer to inviscid 
and viscous solutions. According to this method, we force the quantities 
with subscript v to tend to zero as p + «. Substituting Eqs. (48) and (49) 
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Into Eqs. (44) and (<«^), and letting 6-^0 while holding R fixed, 
we obtain Eqs. (44) and (47). Substituting Eqs. (48) and (49) Into 
Eqs. (44) and (45), taking the limit 6 -► 0 with p fixed, subtracting 
the quantities pertaining to the Inviscid solution, and using Eqs. 
(44) and (47), we obtain 


3% 

+ 1- 

< 

< 

> 

H 

< a 

A 

Rr 


A 

3Z 





Rc 3(1. 




(50) 

(51) 

(52) 

(53) 


The solution of Eq. (51) tending to zero as p -* ® Is = 0. 

It follows that the solutions of Eqs. (52) and (53) tending to zero as 
p -V 00 are 


Vy a^(()>,z) e 


fT p 


(54) 


Wy = a2((j>,z) ^ 

where the real part of VT Is negative, 
into Eq. (50), we obtain 


(55) 

Inserting solutions (54) and (55) 
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3p Rj. 3(j) n 


(56) 


Hence, 


' (’ 


3^1 . VT ; 


/r 1 s* iz 


) e 


(57) 


Writing Eq. (42b) in terms of dimensionless variables, and 
substituting for and q^ from Eq. (48), we obtain 


3U. 3u 3w. ^ 3w 

3z 3Z 3R 3p 


(58) 


Inserting Eqs. (55) and (57) into Eq. (58), we conclude that 


&2 “ O(^) 


(59) 


Similarly, boundary condition (42c) leads to 


“ 0(6) 


(60) 


In view of Eqs. (59) and (60), boundary conditions (41) and (42a) 
become 


A A p 

XF^ + 2u^ = 0(6 ) 

A A 

on R = Rj. 

(61) 

= 0(«2) 

A A 

on R = R^ 

(62) 
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Hence, to order 6, the boi,ndary conditions at the free surface are the 
same as for the inviscid liquid and no dissipation takes place in the 
boundary layer next to the free surface. 

5.3 The boundary layer next to the wetted surface 

We seik a solution valid everywhere except near the free surface. 
To this end, we use spherical coordinates and write the solution in 
the form 


qi “ q^®r Va 6u^(p,a,6)Je^ 


+ [v^(r,a,3) + Vy[p,a,e)]e^ + [w^(r,o,p) + w^(p,a,6)]eg (63) 

where the viscous terms, namely, the components with subscript v, tend to 
zero as p -► <», in which p » (1 - r)/6. The inviscid components, namely, 
the components with subscript i, are solutions of Eq. (47). 

To determine the viscous components, we substitute Eq. (63) into 
Eqs. (44) and (45), subtract the inviscid components, let 6 -» 0 with 

A 

p fixed, and obtain 


1 i 1 

— ^ + — (v sin a) + * 0 (64) 

3p sin OL 3a sin a 33 


A 



AVy - 2fiWy cos a 


^ 2 

— ^ 

3a 3p^ 


(65a) 

(65b) 
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1 


.s 1 ^Py . ^ \ 


XK. 2nv.. cos nt ■ - 


W » • 

T ■ " ^y 

sin a 38 3p 


The solution of Eq. (65a) that tends to zero as p ® is 


Py • 0 


(66) 


from which it follows that the solution of Eqs. (65b) and (65c) is 


SiP SoP 

» c^(a.6)e + C2(a,6)e 


(67a) 


SiP 


S2P 


Wy • i[c^(a,6)e - C2(a,8)e ] 


(67b) 


where 


» X + i2n cos ot 


( 68 ) 


in which i » /-I . Note that the real parts of s-j and $2 must be 
negative for and w^ to tend to zero as p -► ®. Introducing Eqs. (67) 
into Eq. (64), and solving for u^, we obtain 


1 C-. 9S-1 1 i '’'•I. 

“v “ iy s^ ^ sin a 3a ^*^1 * sin a 


i ^"-l. *lP 




1 


S'! 3a 


sin a 38 


Co 3So ^ SoP 

+ — p e 


S 2 3a 


(69) 
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Boundary conditions (40) demand that ^ on r ■ 1 . so that 


■l " 


'1 " \ 


'1 " ^ 


■ 0 

on 

r - 1 

(70a) 

» 0 

on 

A 

r ■ 1 

(70b) 

• 0 

on 

A 

r - 1 

(70c) 


Introducing Eqs. (67) and (69) Into Eqs. (70), and recalling that 
p ■ 0 when r ■ 1 , we have 


+ 6(fl 

+ f2) = 0 

or. 

r - 1 

A 

(71a) 

+ c^ + 

C 2 “ 0 

on 

r ■ 1 

(71b) 

+ i(c^ 

- C 2 ) - 0 

on 

A 

r » 1 

(71c) 


where for convenience we adopted the notation 


c^ 3s. 


3c, 


Equations (71b) and (71c) can be solved for c^ and C 2 » with the result 

Cl * - ^ (vi - 1wi), ^2 “ *^1 “ “ 1 ^''l ^*^1^ 

so that C 2 1s the complex conjugate of c-j. Moreover, from Eq. (68), 

It Is easy to verify that 
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1 

1 1 
1 1 


as 

da 


1 in sin a ^^2 . in sin a 


da 


(74) 


Inserting Eqs. (72) - (74) into Eq. (71a), we obtain 


Ui + ^ in sin a [(^ - ^^'^i “ ^ " f “ * 


1 . 1 


Si S2 


Si S2 


■ ^(i7 - - 1 t<i7 * 




9Wj 


’1 ’2 


's^ S 2 3B 


It will prove convenient to Introduce the notation 


s^ = s e 


-iY 


JY 


(76) 


where 


,-2 . -_2 2 xl/2 1 2n COS a 

s ■ (X + 4n cos a) » Y ■ ^ t^n ^ 


(77) 


Then, the boundary condition at the wetted surface reduces to 


u. . 5 (v^ sin 3y - w^ cos 3y) - 6 (v^ cos y - sin y) 


cot a 


1 dW. 1 ovj o "4 

- « T cos Y * tT- si" y) ♦ « sin Y - 55- cos y) ■ 0 


dw 


3v, 


3W. 


s '3a 


3a 


s sin a 3S 


36 


on 


r - 1 


(78) 



• The Eigenvalue Problem 
6.1 The solution of the equations of motion 

Equations (44) and (47) must be solved simultaneously for the 

A A 

components of the velocity vector and pressure p^ In terms of the 

A 

components of the rotation vector 0. An examination of Eqs. (44) and 
(47), however, reveals that they are not In a form that permits a 
convenient solution, because Eq. (47) Is a vector equation Involving 

A A 

both q^ and p^ , whereas Eq. (44) Is a scalar equation Involving only 

A 

q^. If Eqs. (44) and (47) can be transformed Into a vector equatlo.' 

A A A 

Involving only q^ and c scalar equation Involving q^ and p^ , then the 

A 

vector equation can be solved Independently for q^ , thereupon the 
scalar equation yields p^ . Indeed, taking the curl and divergence of 
Eq. (47), and considering Eq. (44), we obtain the vector equation 


XV X q, - 2n 

3z 


AO A 

+ 2X 6 » 0 


and the scalar equation 


(79) 


A 



respectively. 

Due to the nature of the problem, the use of cylindrical co- 
ordinates is indicated in the solution of Eq. (79). Inserting q.| » 

^ V«l> ® * ®z®z 
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t 


A A 


0p ■ 9^ COS ♦ + 02 ♦ 


0^ ■ - sin ♦ + 02 cos «t» 


(81) 


A A 

e ■ 0, 

z 3 


and equating the coefficients of e^^, e^, and e^ to zero, we obtain 
the three scalar equations 


a(x — - - 2n + 2 a^0r - 0 

R 0 ()) 3Z 3z 


(82a) 


A 3qp 3^7, '' oHi aoa 

A(-t 4 - ^) - 2n ^ + 2A^0, - 


A 3q, 


(82b) 


3Z 3R 


3z 


A 3q. q. 1 3qn a 3q aaa 

A(-^ + - X — ^) - 2 n -^ + 2A^0, - 0 

3z 


3R R R 3(1) 


(82c) 


Equations (82a) and (82b) can be solved slmultareously for 3qp/3z and 


3q,/3z, with the result 
0 


3qr) ■V A n A OM_ A AA AA 

-A^ = AA ^ AA [2n X — ^ + A ^ + 2A(2fl0A - A0. )] 


3 q, A 3q, 


'9 ^9 

3Z + A 


(83a) 


R 3(j) 


3R 


A A A A A 


3q. r A 1 3q, a 3q_ 

[A X — ^ - 2fi + 2A(A0p + 2n0. )] 
3Z + A R 34) 3R ^ 


(83b) 


Inserting Eqs. (83) Into Eq. (82c) differentiated with respect to z, and 
recognizing that 30p/3(j) * 0^, 30^/3(j) ' - 0p, we obtain 
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1 




I 


OA A AA AA 

?2 , aq , aS a q 

A / ^2 . 1 . 1 ^Z\ . _ A 

XJ' ' ■ y-h l "* y~ T •*: -X— ♦ irK ■ - y ) T ’• xy ■ 0 

4fi‘ + A aR R aR R a* 


az‘ 


( 84 ) 


which is a partial differential equation for the velocity component q^ 
alone. Its solution can be obtained by the method of separation of 
variables. To this end, let 




(85) 


so that, inserting Eq.(85) into (84), and dividing through by q^, we 
obtain 




1 . 1 ‘‘^R, . 1 1 ” 'z 

+ X -X-) + xy 9*-J + aJ- 


<|2f_ 


“ a 5 xy [t— (— A-y 

4n^ + A^ ’^R dR^ 


R dR 




T 


♦ 


f^d? 


( 86 ) 


Next let 


42 ♦ a 2 1 


(87) 


where k is a number, and denote 


.2 . .2 

K * J 

+ A'^ 


so that Eq. (87) can be written in the form 


2 

- rf , 0 

dz^ ^ 


( 88 ) 


having the solution 
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(89) 


» 

I 


1 


■ A cosh jz + B sinh Jz 
In view of E!q. (87), Eq. (86) reduces to 


r2 . 1 




^ ft d6 ^ d4)^ 


(90) 


leading to the two equations 


d2f 


— ^ • 0 


(91) 


and 


d^f 1 df» . A 

- ^)fR • 0 


^ R dR 


(92) 


where f Is a number. The solution of Eq. (91) Is simply 


f^ ■ C cos t(|) + D s'l £(J» (93) 

On the other hand, Eq. (92) Is recognized as a Bessel equation having 
the solution 

fp - E J^(kR) + F Y^(kR) (94) 

where and are Bessel functions of order t of the first and 
second kind, respectively. Combining Eqs. (89), (93), and (94), we can 
write 

q^ = (A cosh jz + B sinh jz)(C cos £ 4 » + D sin f(j))[E j£(kR) + F Y^^^R)] 
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(95) 


I 


which represents only the homo«mneous pap t of the solution. Before 
obtaining the particular solution, let us derive the homogeneous parts 
of qj^ and q^. Introducing Eq. (95) into Eqs. (83), ignoring the 
terms in Oj^ and 0^, and integrating with respect to z, we obtain 

^ A 

9d “ “ (A sinh jz + B cosh jz) (-C sin + D cos t^) x 

j R 

[E J^(kR) + F Y^(kR)] + Xk(C cos 44> + D sin f4.)[E J^(kR) + F Y^(kR)]} 

(96) 
and 

q. * - (A sinh jz • B cosh jz){ ^ (-C sin M + D cos l<p) x 

4fr + x'^ j R 

[E J^(kR) + F Y^(kR)] - 2fik(C cos l(P + D sin f(J>)[E J^(kR) + F Y^(kR)]} 

(97) 

A 

where and Y^ designate derivatives of and Y^ with respect to kR. 

Next let us turn our attention to the particular solution, and 
assume the solution in the forr\ 


C,((}))Z 


+ C2(4>)R, 


= C3(4))z + c^((j))R, 


Cg((|))z 


+ Cg(4))R 


(98) 

Introducing solution (98) into Eqs. (82), as well as into the continuity 
equation, Eq. (44), we obtain the set of equations 
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« 


/V 1 '•'Cr . dc^ ^ ^ 

^dr ' " 3T “ ° 

K 


X(c^ - Cg) - 2fiC3 + Z\ 6^ » 0 


1 dc^ dc. 


2 : 


He, + c, f + c. - ^ — i- z - — ^) - zncr + 2X‘0, = 0 


'4 ' '"? ^ ■ '"4 

^ R ^ R d4) 


d(|) 


(99) 


It can be verified that the solution of Eqs. (99) can be written as 


'l = - % • '2 • 0 ■ 's’ «R 


/S A ^ ^ 

C4 = - xe^ , cg = 0 , Cg = X0^ - 2fi0f 


( 100 ) 


A ^ 


Hence, the complete solutions for q^, and q^ are 


Rd ~ ''O ''O 

^ + X^ j 


^ 7 ;^ - (A sinh jz + B cosh jz) (^^ (-C sin + 0 cos ^<1>) x 


[E J^(kR) + F Y^dcR)] + Xk(C cos H + D sin U)H J^(kR) + F Y^(kR)]) 


X0.Z 

4> 


(101) 


Ra = j - (A sinh jz + B cosh jz)( ^ (-C sin + 0 cos H) y 

^ j R 

[E Jfl(kR) + F Y»(kR)] - 2«k(C cos U + D sin J^(kR) + F Y^(kR)]l 


A A /\ A A 

+ x(0pz - e^R) 


(102) 
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q » (A cosh jz + B sinh jz)(C cos + D sin l<ti) x 


[E J^(kR) + F Y^Ckfi)] + (X0^ - 2n0p)R 


003) 


It remains to obtain the solution for the pressure p-| . Introducing 
solution (103) into Eq. (80), and recognizing that F' * e^, we obtain 
the equation 


V^p, “ j (A sinh jz + B cosh jz) (C cos £({) + 0 sin t^i) ^ 

\ 

^ A A 

[E J^(kR) + F Y^(kR)3 - 2nX0^ 

The solution of Eq. (104) can be written in three parts, namely. 

Pi = P2 (R.<I>.z) + P3(R.<^.z) + P4(R) 
where P2 satisfies the equation 
^2 

V^P2 ° J jz + B cosh jz)(C cos £41 + D sin l<i>) ^ 


(104) 


IE J^(k6) + F Y^(kR)] 


A 

p^ is the solution of Laplace's equation, and p^ satisfies 


j'' 

. 1 ^ 

^ P4 “ ''2 ^ ^ ^ 

^ dR^ R dR 


- 2fiX0. 


(105) 


(106) 


(107) 


Writing the solution of Eq. (106) in the form 

= (A sinh jz + B cosh jz)(C cos £cp + 0 sin £(j>)[E* J£(kR) + F* Y^(kR)] 

(108) 
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and recognizing that 


34 »^ 


^ P 2 .2^ 

“ ^ Pp » ''2 * J ?2 

^ 3z^ ^ 


Eq. (106) reduces to 


+ (j^ - 4>] [E* J.(kS) + F* Y (kfi)] 
dR^ R dR R^ ^ ^ 




j [E J^(kR) + F Y^(kR)] 


(109) 


But E* J^(kR) + F* Y^(kR) is a solution of Bessel's equation, Eq. (92), 
from which It follows that 


^ - y ] [E* J,!kR) + F* V,(kR)] 


dR*- R dR 


(j^ - k^)[E* J^(kR) + F* Y^(kR)] = ^ j[E J^(kR) + F Y^(kR)J 

X 


4n ,r 


(110) 


2 2 

Recalling the relation between j and k , we conclude from Eq. (110) that 


E* = 


^ ^ i 


^ ? ? 
X - k*^ 


E = - 4 E, F* = - 4 F 

J J 


Hence, solution (108) becomes 


A 

Pp = - 4 (A sinh jz + B cosh jz)(C cos £41 + 0 sin £(j>)[E J^(k 6 ) + F Y^(k 6 )] 

^ J 
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(111) 


I 

i 

1 


The solution of Laplace's equation can be obtained from the solution 
of Eq. (86) by letting k « J. It follows Immediately that 

Pj ' (A-| cosh jz + sinh jz) (C^ cos £<j> + sin f(|») x 

[E J^(jR) + F Y^(jR)] (112) 

On the other hand, the solution of Eq. (107) can be written as 

P4 ^ nAe^R + G (113) 

where the meaning of the constant G will become evident later. 

Combining Eqs. (Ill), (112), and (113), we obtain the pressure 

A 

Pi “ " j (A sinh jz + B cosh jz) (C cos £(}> + D sin t<t>) x 

[F J^(kR) + F Y^(kR)J + (A^ cosh jz + B^ sinh jz) x 

(C^ cos + 0^ sin i<t,) [E^ J^(jR) + F^ Y^(jR)J - | nAB^R^ + G 

(114) 

Upon Introducing solutions (101) -(103) and (114) Into the moment 
equation, Eq. (46), and performing the appropriate Integrations, It is 
possible to eliminate the constants O2* and from these solutions, 

AAA A A A 

where we note that G-j , 02* and are related to 0^, 0^, and 0^ by 
Eqs. (81). To this end. It is more convenient to work with the body axes 
X, y, z. The Inertia dyadic T has the principal diagonal elements 

AAA 

111, I22. I33 and the off-diagonal elements - I-|2i - 1^31 - I23’ 
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" ^21 " ” H 2 * ’ ^31 “ " hs’ " ^32 " “ ^23* the body axes 

to coincide with the principal axes, the products of Inertia reduce to 
zero. Consistent with previous assumptions, we shall Ignore the term 

A A 

Rg In Eq. (46). In view of the solutions for and p^ , Eq. (46) can be 

A 

written In the form of three simultaneous algebraic equations In 0^, 
e^, and 0^. These equations can be conveniently displayed by Introducing 
the notation 


M 


11 


S A^z^(cos^(j) - s1n^(j>) + (n^ - X^)6^ s1n^4> - 2fiX(R^ - 2z^) > 


m. 


s 1 n 4 > coS ( j ) Jdm . 


M 


12 


,.0/sA ^ 0^0 aaaA P aaaO ^ 

^[(fi R - X R - 2X z ) Sln 4 ) cos(j) + 2ftXR s1n^(j) - 2nXz ]dm^ 


ni / 


f AAA A A A 

Mi 3 = I XRz (X cos4> - n s1n4>)dm^ 


= f [2«XR^ cos^(t> - 2lUz^ - (n^ - X^)6^ s1n4> cos4)]dm^ (115) 

'' m » 


'22 


^ [_X^z^+ (n^ - X^)R^ cos^(j) + 2nxR^ s1n(() cos^'Jdm^ 


m. 


^23 ^ j'' ^ ^ coS(J))fiz dm^ 


M 


31 


[- is? - X^)coS(|> + Zrs s1n(|)]Rz dm^ 


A A _ AA . A 


-'m. 
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A ) sirv)) - 2nA cosi}>jRz dm 


”33 “ 


aO aO aO a 

(fl - A ) R dm, 




where (1,j = 1, 2, 3) can be Identified as terms involving the 

^ A 

nonhoinogeneous parts of and p.| , and 


A /V 


{( + Azq. + 2flzqo) COS(j) + [z — R -;r- - RR,) 

AC. m K rv _ '' 


m^ d(p 


4> --R^ 

- 2nzq^3 sin<{)}^ dm^ 


^P^ ys ^P^ A .AA AA 

3^ 3Z 


H 


jjh 

Jm. R 


A A 

Z 


A AA AA 


9Pl AAA AAA A 3P] . 3P. - 

o - i^iy^ — - * Xzq. + 2nzqp) sin<t) - [z — - R — + A(zqp - Rq ) 
^ 'm, R 8* ♦ ^ 3R 3z 


- 2fizq^] coS())}^ dm^ 


(116) 


m 


3Pi AAA AAA A 

<3T * * ^“'"'r’h “"’e 


represent terms involving the homogeneous parts of q^ and p^ , where 
the subscript H indicates that only the homogeneous parts of the solutions 

A ^ 

for "q^ and p^ are to be included in the integrands. Considering the above, 
Eq. (46) can be ‘ Uten in the form of the three scalar equations 


31 


^AA A A AAA AA A A ^ ^ 

+ xni^2 ■*■ ^ hz ^ ^”^^3 ' ^22^ ^ ^2^®2 

a9a AAA A A 

* (- » 1,3 ♦ 2Anl23 + ", 3)03 - H, 

aOa a a a a a a aOa AAA A A 

[“ X 1^2 '*’ XnCl'ji •■ ^23^ ^ ^22 ” Xfll^2 ^22^®2 

(117) 

A O A AAA A A A 

+ (- X‘l23 - 2XfiI^3 + M23)03 ■ H 2 



Xftl23 + 


«31>®1 


A O A 

+ (- XM 


23 


+ Xnl^3 + M32)02 


(x^i 


33 


^33^®3 


A 



A 

Equations (117) can be solved for 0 ^ (1=1, 2, 3) to eliminate these angles 

— A 

from and . 


6.2 Satisfaction of the boundary conditions 

A 

The solutions q-j and contain 13 constants of integration in 
addition to the numbers j and t and the eigenvalue X. To determine 
these quantities, we must invoke the boundary conditions. The kinematic 
boundary condition, Eq. (61), merely yields the function provided 
the solution is known. On the other hand, boundary conditions (62) and 
(78) can be used to determine the quantities in question. 

At the free surface boundary condition (62) reduces to 

p^ = 0 on R = (118) 
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( 


Fiom Eq, (114), we conclude that boundary condition (118) Is satisfied 
for all ip and z if 


J.(kL) 

' E . F, E, . G - i sne R‘ 

I V / jfi \ I c 2 C 


(119) 


Y^(k6c) ' Y^(JR^) 

From the structure of the solution (114), however, we observe that the 
constants E and E^ can be absorbed by A and B and and 6^, respectively. 
To this end, and for convenience, we let arbitrarily 


E 




1 


so that the pressure becomes 


(120) 


Pi ■ - j (A slnh jz + B cosh jz) (C cos £(|) + D sin i(p)>t 

[Y^(kRc) J^(kR) - J^(kR^) Y^(kR)] + (A^ cosh jz + slnh j$) x 
(C^ cos U + sin tp) [Y^(jR^.) J^(jR) - J^(jRc) Y^(JR)] 

1 aO aO 

+ ^nX02(R^ - R^) (121) 

Moreover, the velocity components become 

% “ - (A slnh jz + B cosh jz){ (-C sin £(j> + D cos ip) x 

^ ^ i R 

[Y^(kR^) J^(kR) - J^(kRj.) Y^(kR)] + U(C cos tp * D sin f<^) x 
[Y^(kR^) J'(kR) - J^(kR^) Y'(kR)] - Xy) (122) 
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• I 

( 


— 77 “ jz + B cosh jz) {^ (-C sin «,(j) + D cos 4(|))x 

* j R 

[Y^(kR^.)J^(kR) - J^(kR^)Y^(kR)] - 2J2k(C cos £<» + D sin £o)x 

q^ » (A cosh jz + B sinh jz)(C cos £4» + D sin £^)[Yj^(kR^)J^(kR) 

The boundary condition at the wetted surface Is as given by Eq. (78). 
Satisfaction of this condition at every point of the boundary should yield 
the effects of viscosity on the stability. 

7. Semi -Analytical Solution 

Because the solution (121) - (124) Is In terms of cylindrical coordinates 
with the origin at the satellite center and the boundary condition (78) Is In 
terms of spherical coordinates with the origin at the center of the spherical 
tank, no exact solution Is possible by separating variables. In view of 
this fact a semi -analytical numerical soltulon has been attempted. 

The semi-analytical solution consists of solving the problem In two 
stages. In the first stage a solution of the eigenvalue problem is obtained 
by regarding the liquid as being entirely Inviscid, and In the second stage 
the Inviscid solution is perturbed, so that the perturbed motion satisfies 
the viscous boundary condition at the wetted surface. The Inviscid solution 
was attempted by a Galerkin procedure, whereby a solution was assuned In the 
form of series of admissible functions as follows: 
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( 


• y.^ A^Qpj.(R,4.,z) , 


q. ■ y c Q (H,(j),z) , 

2 S-1 * “ 


A 



n 

y B 0. (R.^.z) 
s-1 * 

n 

y D P (R.^.z) 
s-1 ^ * 


(125) 


where Qn,t Q . , . Q„. and P, are functions similar to solutions (121)-(124), 

KS 9$ ZS S 

but satisfying the boundary conditions at both the free and wetted surfaces. 

We note that these functions no longer satisfy the continuity equation and 
the equations of motion. 

Introducting tqs. (125) Into Eqs. (44), (46), and (47), multiplying 
Eqs. (44) and the three components of Eq. (47) by ^Rt’ %t* ^zt and P^, 
and Integrating over the volume of liquid, we obtained an algebraic eigenvalue 
problem of order 4n+6 In the coefficients A^, B^, C^, D^, the angles of 

^ A ^ A 

0^ , 02 and with X as a parameter. The eigenvalue problem and the 

A A ^ ^ A ^ A A A 

angular velocities 0)2 * Xuij^e^, and ■ ^1^63 can be written 

In the form 


X EX + Fx - 0 (126) 

A solution of the eigenvalue problem (126) was attempted by the QR 
method but, in spite of considerable time and effort, the results proved 

A 

unsatisfactory. In particular, the obtained eigenvalues X have real parts, 
which is in contradiction with the expectation of a purely oscillatory 
motion for an Inviscid stable problem. The conclusion that can be drawn 
from this is that a Galerkin approach is not feasible for this problem be- 
cause of the difficulty in selecting a satisfactory set of admissible functions. 

Because the semi -analytical approach described above did not lead to 
satisfactory results, a completely numerical solution of the Inviscid pro- 
blem with the modified boundary condittlons may be attempted by means of fi- 
nite differences or finite elements. Although that approach may permit 
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a ready extension to other tank goometrles. the work promises to be 
extremely tedious, because a three-dimensional finite-difference grid 
or a three-dimensional finite element must be used. We note. In passing, 
that truly three-dimensional numerical solutions of fluid dynamical pro- 
blems are scarce. 
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